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A  finite-element  methodology  has  been  presented  to  solve  the  one -dimensional 
couple  dynamical  thermoelasticity.  The  spatial  variable  is  discretized  by  the 
Galerkin  scheme  and  the  central  difference  scheme.  The  time  variable  is  dis¬ 
cretized  by  the  Galerkin  and  the  backward  difference  scheme.  It  was  found  that 
the  two  spatial  discretization  approaches  give  close  comparison  in  accuracy  whei 
used  in  conjunction  with  the  backward  difference  in  time.  To  simulate  the  wave¬ 
like  behavior  in  micro-time  scale  becomes  expensive  but  there  is  no  technical 
difficulty.  To  simulate  the  macro-time  phenomena,  the  micro-time  behavior  woul< 
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Two  test  problems  are  presented  to  show  the  micro-  and  the  macro-time  behaviors 
of  the  thermoelastic  slab.  The  responses  to  a  realistic  set  of  pressure  and 
temperature  pulse  inputs  simulating  the  interior  ballistic  type  of  action  are 
also  presented. 
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A  FINITE  ELEMENT  SOLUTION  OF  THE  COUPLED  DYNAMICAL  THERMOELASTICITY 


ixtrotvctitn* 

The  coupled  themce  Lascicity  p rcblam  has  seen 
-  reared  by  many  investigators  in  she  past.  A  partial 
listing  of  the  relevant  capers  are  given  in  the  refer¬ 
ence  '1-53#  most  of  these  papers  leal  with  the  problem 
of  a  seni -infinite  space  and  a  few  contain  a  comprehen¬ 
sive  numerical  treatment.  In  many  applied  problems  the 
coupling  cam  in  the  energy  balance  equation  is  ignor- 
able  cue  to  the  snailr.ess  : :  its  magnitude.  Ir.  other 
problems  such  as  that  tf  the  thermal  stresses  in  a  gun 
barrel,  cue  to  the  hluh  rats  o :  loading  the  coupling 
.oust  be  Jcect  until  the  results  of  the  analysis  shew 3 
otherwise.  Also,  the  martial  effects  might  be  impor¬ 
tant  and  shot  Id  be  dept  ir.  t.te  analysis  for  such 
problems.  The  particular  mathematical  model  to  be 
dealt  in  this  paper  is  a  boundary  value  problem  for  a 
pair  of  partial  cifferential  equations  governing  the 
stress  and  the  temperature  fields.  The  boundary  con¬ 
ditions  correspond  to  tie  time  dependent  stress  and 
temperature  acting  on  one  surface  and  the  free  traction 
and  ambient  temperature  on  me  other  surface. 

A  finite  element  approach  is  used  in  discretizing 
the  spatial  as  well  as  the  time  variable.  In  each 
discretization  a  choice  of  two  different  approximation 
is  provided  to  furnish  a  comparison.  Although  the 
problem  analyzed  is  or.e-dir,ensior.al  in  space  the  method¬ 
ology  developed  can  be  transferred  to  the  problem  of  an 
annul us  without  too  much  chance  in  its  principles  and 
cocr.ni cues .  Ir.  the  next  section  the  mathematical  medal 
will  be  giver,  in  a  general  form  first  and  then  special¬ 
ized  cc  a  slab. 


where  •  is  the  Laplacian ,  •.  ar.c  u  are  me  _ane  s  con¬ 
stants,  £  is  the  body  force  vector,  c  is  the  density 
and  v  is  a  constant  given  by 

y  *  (3A-*-2u)a  12) 


where  a  is  the  coefficient  of  thermal  expansion.  The 
last  term  in  the  left-hand  side  of  Zq.  Cl)  is  the  force 
per  unit  volume  cue  to  the  temperature  gradient  in  the 
medium.  In  an  ur.ee uo led  oroblem  the  tencerature  field 
is  given. 

The  second  aquation  of  the  pair  of  partial  dif¬ 
ferential  equations  of  the  moadel  is  the  energy  balance 
as  given  by  ill 
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where  <  is  the  dif fusivity ,  1  is  the  heat  source  ar.c 

n  =*  vT  /cc<,  T  is  the  ambient  temoerature. 

00 

Assume  that  the  body  force  £  and  the  heat  source  2 
can  be  dropped  and  that  the  displacement  vector  u  is 
of  the  form 


u  «  grade  +  curl  0  (4) 

Substituting  (4)  into  (1)  and  ; 3 )  yields  a  pair  of 
•couple  PD  £  in  $  and  T  and  an  uncoupled  ?D2  in  .  des¬ 
cribing  the  shear  wave  propagation.  These  equations 
are: 


(5) 

(6) 


(T) 


where  m  *  y  ( i-*-2u )  ,  Cf  c  .X*2u)/o  and  Zi,  *  u/c  ,  2 .  and 

1  *  j. 

and  C2  are  the  diiatational  and  the  shear  wave  speed 

respectively.  Since  the  temperature  fields  has  no 
interaction  with  the  shear  wave,  the  subsequent 
analysis  will  direct  its  attention  to  (5)  and  (6). 

To  specialize  the  pair  of  equations  to  a  slab  we 
can  write  (6)  in  the  form 


and  (5)  as  follows 


(S) 


••ttTKSXATrCAC.  MCOEL 

let  the  temperature  in  the  medium  be  T  and  its  dis¬ 
placement  be  fenoted  by  u,  then  the  equation  of  motion 
of  the  me ^ 1 urn  is  given  in  a  vectorial  form  as  follows. 

^T*u*  '  *■  +'j)  grad  div  u  *  £  -  ygradT  *  oU  (1) 
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Since  it  is  more  convenient  to  use  the  stress  instead 
the  potential  :  as  the  variable  in  one-dimensional 
problem,  we  introduce  the  constitutive  equation  for 
a  thermoelastic  medium. 
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THE  SPACE  COMA IN 


where  c  is  the  uniaxial  stress. 

Jsir.g  f  1 3 )  to  -aiimir.ata  5  from  '3)  yields 


To  discretize  the  spatial  variable  in  (15)  and 
(17)  we  used  two  approaches ,  i.e.  the  3a lark in  proce 
dure  and  the  central  difference  scheme. 


The  Galerkir.  Procedure 


In  this  procedure  we  snail  assume  that  both 
approximate  solutions  of  (16)  and  (17)  *  and  I  car.  be 
expanded  in  terms  the  same  spatial  approximate  function 
v.  (x)  in  the  following  way. 


The  boundary  conditions  are  riven  as  fellows, 
x  =  o,  c  (2,t  >-f*  (c)  *  0 

7  'T.'0,t)-T  )-(g*rC)-T  )  *0  (13) 


5  *  /  v .  (x)c .  (t  x  a  (t) 

>.11  co 
1*1 


)  v. (x)T. (t) 
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1*0 


where  i  is  t re  thickness  of  the  slab,  6*  *  h^/!<, 

6*  =  h?/K,  h.  and  h,  are  the  heat  transfer  coeffic 
and  the  conduct ivi tv  of  the  material. 


where  >3=r.+l ,  n  is  the  number  of  internal  nodes  and 
P.  (x)  is  defined  as 


1-1  — ^  i  ,  1  x-x.  *  <  h 


non- j imsns ::nal  Eg cat ions 

To  facilitate  computation  the  following  rvon-di- 
mensior.al  quantities  are  introduced 


After  converting  the  equations  into  the  non-dimer.sional 
font  the  bars  are  dropped  for  ease  of  writing. 

Thus,  after  some  algebra,  we  have 


-  *  0  ,  |X-x..>h  i=0,l,2, — 7*  ( 21 ) 

where  h  is  dimensionless  length  of  the  interval.  A 
graphical  illustration  of  the  function  in  (21)  is  given 
in  Figure  1. 

For  convenience  we  introduce  the  symbols  and 
L2  such  that  the  PDE  in  (15)  and  (17)  can  be  abbreviated 


L,  (o,T)  *  0, 


~  _i  *  =o 
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L2(o,T)  *  0. 

Then  the  residue  and  3.,  can  be  defined  to  be 


2  ’ xx  *  1 '2 '  ' t  ^l7t 

_  „  .  ,  ,  ,3\+2u,  *  , 3X+2j .  „ 

sre  _*C.  «./ < /  \j*  v— — )a.  a2=  (  ^  )*  T^. 


R1  ‘  and 


*2  =  -2w..:. 


For  the  boundary  conditions  we  abbreviate  ‘'13)  ar.d  (13) 


The  boundary  conditions  reduce  to 


Bj  (o)  =  0,  B2  !T)  *  0  and 


x  *  0.  7(0,0  -  f'O  •  0 


S2(o)  -  0,  32(T)  *  0 


Tx  -  i,T  -  3.,  -  0 


X  -  1  7(1,0  *  0 


In  the  Gaierking  procedure  we  require  that  the 
weighted  residue  to  vanish  with  the  weighting  function 
taken  as  v  .  ix'  ,  i  .e . 


T  -  3,T  •  3 
x  2 


where  22'f*i,  f<0«  f*/('-+20,  g- 


/  S.i . (x)dx-l  i  *  1,2 
n  1  1 


i  -  3,1,2, 


For  1  »  1,  (26)  leads  to  after  integration  by  parts  and 
some  algebra. 
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v  .  fx)*y  .  (x)7.  't)dx 

1  2  i  ii>  1  3  1 


+3.  j  >  a. !x)5 (x) g(t)dx 
0  i-O  1 


i  '  x)  i .  'x)  3  .  (c)  dx—  -V,  ^  (x)b.(x)o  (t)dx 

111  l'.3:d 


-3,  ;  .)  a,  (x)3  (x-l)T;  (t)dx 

0  1=0 


i!  tx)  ti,.  (x)T.  !t) dx! 


Performing  the  necessary  integration  in  (27]  and  ,32] 
yields  the  x-discretized  equation  of  (11) . 


To  treat  the  term  in  the  bracket  in  (27)  the  ther¬ 
mal  boundary  renditions  in  '21 ;  and  25)  will  be  used, 
let 


(1+S,h)  -1 


i  -1  2-1 


(3)  5  3,  (T) 


-1  2  -1 


5>:<X)T.  (e)-s.  f  3.  (X)?.  (t)*J,?(t)  (23) 
i-0  *  A  i-0  1  1 


sc;  3  3,(7) 


=  /  j!  !x)T.  (t)+3  y  ii.(x)T.'t)  (29) 

i=0  11  2  iio  1  1 

"*e  demand  chat  the  weighted  residue  of  R(0)  and 
3(1)  vanish.  This  statement,  when  (23)  and  (29)  are 
i3sc,  ran  be  tasted  inti  the  following  form  by  intro¬ 
ducing  a  delta  function  in  each  term. 

1  N 
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-3,  ]  i  l .  (x)  ii .  (x)T.  (t)  3  (x)dx 
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*3.  /  u  'x)5 (x)g(t)dx  -  0 
*  0  3 


u!  (x!  (x)T. (t) 5 (x-lldx 
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*3,  ,  3  ■ x) 5 (x-l)T. (t)dx  « 

‘  :  i-o  1 


■J.  (X)  »  fx)T.  (t)  !  - 
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1 

For  i  »  2,  /  R  '|i.(x)dx  »  0 

0  : 

9y  a  similar  calculation  as  above  we  obtain  the  x-dis- 
cretiaed  equation  of  (12) 


-y  (A21  £  +  [B2J  £  +  >2  [HI  T  T  -  f2  (35) 


where  (A.  1 ,  [A  |,  etc.  are  the  coefficient  matrices  in 
(33)  and" ( 34) . 2 

Central  Difference  Scheme 

In  this  scheme  the  approximating  functions  -j*  (x) 

are  quadratic  instead  of  linear  as  in  the  case  of  the 
Gaierkin  scheme.  The  weighting  function  is  taken  to  be 
•5  (x-x^ )  . 

Let 

i+1 

a(x,t)  -  l  J>..  (x)a  .  (t)  (36) 

i-l  -  3 


(35) 


and  setting 


1 

/  S<0)  <5  ( x-0 ) dx  -  0 
0 

result 


_  (-T_lwri>  -  3j_To  +  3^-0 


from  which  we  can  calculate  T  , 

T-1  *  V2h8lTo  +  2hV 


(42) 


.Substituting  (42)  into  (41)  yields 


1  -281- 

U*h8.)-2T.  1  +  (1+X,  X,)T  -  —■  g-X,  f  (43) 
CJ>2  o  11  1  2  o  Ch  1 


Similarly,  setting  /  R(l) 5 (*-l)dx“0,  we  have 


^j(-2Vl*2«l*«2>‘)^l*(l*X1X2)TH  -  0  (44) 

cn 

From  (40.  (42)  and  (44)  we  finally  obtain  the  x-dis- 
cretited  equations. 
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(45) 


(46) 


equations  (45)  and  (46)  can  be  abbreviated  as 
follows. 


(47) 


JL  r. 

h2(*JJfl7’C'2la*x2lc)  I  -  f2 


where  (E  ] ,  [E  ],  etc.  are  the  respective  coefficient 

matrices lin  (45)  and  (46).  Notice  the  equations  in  C  K 

(4S)  are  of  the  same  format  as  those  in  (35),  only  the  it  +  2  ~^i-l+"-i+l  +“S  6  ^i-l*4^l+^i+l* 

coefficient  matrices  are  different. 


THE  TIME  DOMAIN 


(51) 


Since  Eqs .  (35)  and  (47)  are  of  the  identical 
format,  only  (35)  will  be  discussed  in  detail. 
Equation  (35)  can  be  further  abbreviated  as  follows. 
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To  further  abbreviate,  write  (48)  as 


Mjr  +  Cy_  +  -  £  ■  o 


(49) 


where 


(N+n) , 


other  symbols  are  self-iden¬ 


tifiable. 


In  discretizing  the  time  domain  we  shall  develop 
the  recurrence  schemes  via  i)  Galerkin  and  ii)  back¬ 
ward  difference  schemes. 


Rearranging, 

(M+  +  ^|i)z.+i-(2M-  fdt/Oy.-MM- 


-  ^(F.  +4F.  +F  ) 

6  -3.-1  —1  —1+1 


(52) 


Approximating  the  First  Step 

For  the  recurrence  scheme  to  proceed,  both  ^  11113 
^  must  be  known.  An  approximation  of  £  i3  obtained 

by  linear  interpolation  of  the  function  £(t)  over  the 
interval.,  i.e. 


2.(t)  -  NQ(t)2o  +  ^(t)^ 


(53) 


Aqain  requiring  the  weighted  residue  to  vanish,  i.e. 
At 

/  R  N,  (t)dt  »  0 
0 


we  haver  after  a  lengthy  calculation. 


<7  +  f  K>*i  *  -  T  *  ( 3  ’  l 


(54) 


Backward  Difference 


Galerkin 


Expanind  j;  in  an  infinite  series  [6] 


In  this  scheme  a  quadratic  interpolation  between 
two  successive  time  increments  is  used  and  the  weight¬ 
ing  function  is  J (t-t^+^) . 

The  vanishing  of  the  weighted  residue  yields 


1 


l  yt>  ^ 

i-0 


(M+3-|^  Cj^^-flM+Idt+C^+lM+^+K);^-  F^+1  (55) 


f 

where  N.  (t)  -t 

L 


t-t 

"I* 

1 

at 

- 0 ,  ! 

it-tj 

t-t. 


>  it 


<  it 


Abbreviating  (49)  as  L(jj)  «  0  and  define 
3(£)  I  (.(£).  Getting  the  weighted  residue  vanish 
and  performing  the  necessary  integrations,  we  have 


COMPUTATION  AND  DISCUSSION 

In  the  previous  sections  methodologies  based  on 
finite  element  approach  have  been  developed  to  treat 
the  boundary  value  problem  of  the  coupled  PDE  of 
dynamic  thermoelasticity.  Some  computational  results 
will  be  discussed. 

Different  combinations  of  the  spatial  and  time 
discretization  were  implemented.  In  the  graphs  pre¬ 
sented  at  the  end  of  this  paper  a  notation  such  as 
"CD/BK"  would  mean  central  difference  in  spatial  dis¬ 
cretization  and  backward  difference  in  time  discreti¬ 
zation.  The  notation  GX  stands  for  Galerkin. 

The  parameters  of  the  mathematical  model  are  only 
suggestive  values  which  would  provide  some  idea  about 


a  realistic  problem.  Vi  thin  this  objective,  we  have 
chosen  the  following  set  of  parameters  X,  »  1.0, 

X,  -  10'3  (to  10'2),  T  -  100*F,  3  -  S,  »  .4. 

2  0  12 

There  are  two  scales  of  time  for  the  thermo- 
elastic  problem.  One  time  scale  corresponds  to  the 
time  of  travel  through  the  slab  of  a  dilatational 
wave.  The  other  is  related  to  the  time  of  diffusion 
phenomenon.  The  difference  in  magnitudes  of  these 
two  times  span  five  orders  of  magnitude.  The  wave 
travel  time  is  unity  in  the  non-dimensional  time 
defined  in  Section  3.  The  "diffusion"  time  is  about 

103  non-dimensional  units.  Therefore,  basically  there 
exist  two  types  of  responses  in  this  type  of  thermo¬ 
elastic  problem,  which  we  shall  call  micro-time  and 
macro-time  behavior  respectively.  Thus,  the  time  in¬ 
crements  it  used  in  the  computation  also  differ  by 
similar  orders  of  magnitude. 

Three  test  problems  will  be  discussed:  (1)  re¬ 
sponses  to  a  unit-step  stress  f(t),  (2)  responses  to 

a  unit  temperature  step  g(t),  (3)  responses  to  both 

a  stress  and  a  temperature  pulse,  simulating  the  gas 
pressure  and  temperature  in  a  gun  barrel. 

A  Modified  'Jnit-Step  Stress  Input,  f  (t) -1-e' 

'Jncoupled  Problem.  The  theoretical  response 
(solution  to  the  wave  equation)  of  a  slab  to  a  unit- 
step  stress  input  imposed  at  x  >  G  is  given  in 
Figure  3  by  the  rectangular  waves  for  x  ■  0.24.  This 
theoretical  answer  is  used  for  comparing  the  accura¬ 
cies  of  the  responses  to  the  modified  step  computed 
by  the  several  combinations  of  It  and  lx  as  exhibited 
in  the  tabulation  in  Figure  3.  This  figure  shows 
that  decreasing  it  and  &x  improves  the  accuracy. 

The  asiplitude  of  the  response  using  lx  -  .05  and 
It  -  .002  matches  the  theoretical  solution.  Other 
combinations  show  deterioration  of  accuracy  after  one 
or  two  cycles. 

Figure  4  shows  the  comparison  of  the  quality 
of  the  computational  results  based  on  different 
discretization  schemes.  It  is  3hown  that  CD/BK  and 
GK/BK  give  quite  close  results.  In  other  words,  the 
spatial  discretizations  by  CD  and  GK  do  not  yield 
significantly  different  results.  The  backward  scheme 
in  time  introduces  some  artificial  "damping",  causing 
inaccuracy,  whereas  the  Galerkia  scheme  in  time 
generates  oscillations  of  the  amplitude  about  the 
exact  value.  These  oscillations  can  be  reduced  bv  de¬ 
creasing  both  it  and  ix.  Figure  4  shows  that  using 
lx  »  .05  and  It  •  .02  reduces  the  oscillations  when 
compared  to  larger  ix  and  it.  The  graphs  showing  the 
results  of  larger  increments  are  omitted  for  the  sake 
of  space. 

Coupled  Problem.  Figure  5  shows  the  stress  and 
the  temperature  response  dua  to  a  modified  unit  stress 
applied  at  x  »  0  and  zero  temperature  input  at  x  •  0 
with  coupling  parameters  *  .01  and  X ^  «  .1. 

It  is  shown  the  propagation  of  the  stress  wave  and 
also  the  temperature  wave  due  to  the  coupling  effect. 
The  solid  lines  show  the  wave  front  at  various  time 
instants.  The  dotted  lines  show  the  temperature  wave 
being  synchronous  with  the  stress  wave.  The  magnitude 
of  the  non-dimensional  temperature  is  1/10  of  that  of 
the  stress.  This  ratio  is  clearly  determined  by  the 
coupling  parameter  X  («0.1).  It  is  also  observed  that 
at  the  time  corresponding  to  that  of  the  wave  travel 
through  the  slab  (micro-time) ,  no  diffusion  affact  is 
observed. 


A  Modified  Step  Temperature  Input.  q(t)  »  30(l-e'3l:) 

Uncoupled  Problem.  Figure  8  shows  the  temperature 
profile  for  various  times  due  to  an  input  of  a  modified 
temperature  step  with  amplitude  30,  for  the  uncoupled 
equation  (X^  *  X^  =  0) .  Since  the  thermal  problem  has 
a  long  time  scale,  the  result  demonstrates  that  using 
it  ■  20  and  500  yields  no  significant  difference. 

Figure  7  shows  the  time  response  at  two  different 
field  points.  On  the  same  figure  the  instability  of 
the  Galerkin  time  discretization  is  demonstrated. 

Coupled  Problem.  Figure  3  and  9  are  the  results 
in  responses  for  the  coupling  specified  by  X^  *  .01 

and  two  values  of  X^  «  0.1  and  10.0.  It  is  observed 
that  the  temperature  responses  show  significant 
difference  when  the  coupling  parameter  X^  is  10.0, 
which  is  an  upper  bound  value. 

Dual  Pulses,  f(t)  ■  .1  te  ~3t,  q(t)  »  8.1548  te  '3t 

Simultaneous  inputs  of  gas  pressure  and  tempera¬ 
ture  as  approximated  by  the  given  set  exponential  func¬ 
tions  yield  responses  at  x  »  .2  in  micro-time  as  shown 
in  Figure  10.  It  Is  observed  that  the  wave  front 
arrives  at  x  «  .2  at  t  ■  .2.  The  successive  cycles  of 
wave  travel  are  depicted  by  the  ripples  oscillating 
about  a  mean  curve  following  the  general  input  pulse 
wave  form. 

Figure  11  shows  the  micro-time  response  in  the 
temperature  and  the  stress.  For  X  «  1.0  and  x2  -  .01 

the  top  figure  shows  that  the  coupling  effect  on 
stress  is  minimal,  whereas  the  lower  figure  shows  that 
the  coupling  has  a  strong  effect  on  the  wave-like  tem¬ 
perature  profile. 

Figure  12  and  13  show  the  macro- time  temperature 
vs.  time  and  the  temperature  profile  under  the  action 
of  a  dual  stress  and  temperature  input.  Notice  that 
in  macro-time  the  stress  remaining  is  low.  This  is 
due  to  the  fact  that  the  slab  is  free  to  expand. 

CONCLUSION 

The  theory  and  the  implementation  of  a  finite 
element  methodology  in  solving  the  problem  of  the 
couple  dynamic  thermoelastic  slab  has  been  established. 
Test  problems  of  a  unit-step  stress  input  as  well  as  a 
unit-step  temperature  input  are  used  to  exhibit  the 
two  different  types  of  responses,  namely,  the  micro- 
and  the  macro-time  behavior  of  the  slab.  In  the 
spatial  discretization  both  the  Galerkin  scheme  and 
the  central  difference  are  satisfactory  whereas  the 
backward  difference  scheme  is  preferred  in  the  time 
discretization.  A  realistic  set  of  the  stress  and  the 
temperature  inputs,  simulating  the  interior  ballistics 
of  a  gun  barrel,  is  used  to  generate  responses. 
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Figure  5.  -  Stress  and  Temperature  vs  x 
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